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Abstract. We propose a new method for implementing process tomography 
that is based on the information extracted from temporal correlations between 
observables, rather than on state preparation and state tomography. As such, the 
approach is applicable to systems that are in a mixed state, and in particular to 
thermal states. We illustrate the method for an arbitrary evolution described by 
Kraus operators, as well as for simpler cases such as a general Gaussian channels, 
■ ^ ^ ■ and qubit dynamics. 
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1. Introduction 

Quantum process tomography deals with estimating the dynamics of unknown 
systems. It has long been studied due to its fundamental importance in the fields 
of quantum communication and quantum computation. 

The common approach for implementing quantum process tomography is based 
on applying the dynamics on each element of a complete set of input states, and 
then performing tomographic measurement of the output states. This allows one to 
completely reconstruct the superoperator representing the dynamics [T|. In recent 
years several improvements have been proposed for this method [21 [SI H] i which reduce 
its complexity in some cases. 

Another approach is based on applying the dynamics on random states and 
comparing the results with the theoretical output state that would have been received 
had the transformation been purely unitary [5]. The major advantage of this method is 
that it is efficient (in the sense that it scales polynomially with the number of particles 
in the system). However, it does not reconstruct the superoperator representing the 
channel, but rather estimates the strength of the noise. 

In the present work we propose a different method for preforming process 
tomography of an unknown evolution, which is based on the information extracted 
from temporal correlations between observables, and as such does not require state 
preparation and state tomography. In fact, it can be used with an unknown mixed 
states, and in particular with thermal states. As such it is potentially useful in "hot" 
systems whose state cannot be controlled with current technology, or cooled to the 
ground state. 

The proposed method is based on two key points: (a) we utilize a proper set of 
temporal correlations between observables at two consecutive instances of time, say 
t 2 and t\, ti > ti, that encode the dynamical evolution of the system, and (b) we 
employ weak measurements, rather than ordinary disturbing measurements, in order 
to measure such correlations. We cannot use ordinary measurements to observe the 
above temporal correlations, because the interaction with the system at the earlier time 
will disrupt the correlations. Nevertheless as we show, by using weak measurements 
which barely affect the system such correlations can be measured. 

Our method is applicable for general discrete systems whose evolution is described 
by a set of Kraus operators, as well as continuous systems whose evolution is linear 
with respect to a set of operators. The rest of the paper is organized as follows: In 
the next section we define a more general temporal covariance matrix which will be 
sufficient for determining the above type of general evolution. We then show how the 
temporal covariance can be measured using weak measurements. In section 3 we show 
how to find the Kraus operators for general discrete systems, in section 4 we show 
how to find the dynamics for Gaussian channels and in section 5 we demonstrate the 
proposed method on qubit dynamics. 

2. Measurement of a temporal covariance matrix 

Let us begin by defining the two-point temporal covariance matrix: 

a t] {t u t 2 ) = ({BiitJ-iBiihV^jih)-^ (ta)>}> 
= 2Re[(B i (t 1 ) B 3 (t 2 )) - (Bifa)) (B 3 (t 2 ))] , 



(1) 
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where {Bi} are Hermitian operators. This unusual definition of the covariance matrix, 
where the measurements are taken at different times, proves to be essential for our 
method. 

In order to measure the covariance matrix one has to measure the observables 
Bi (t) and the correlations {Bi (ti) , Bj (£2)}- While the measurement of the 
observables is straightforward, the measurement of the correlations is not trivial and 
the rest of this section is devoted to that. 

The general scheme for the measurement of a single correlation {Bi (t\) , Bj (£2)} 
involves two measurement devices: one measures Bi at time t\, and the other measures 
Bj at time t 2 . In order to measure the correlation one has to measure a joint operator 
of both measurement devices. Our method can be implemented in numerous ways. 
For simplicity we elaborate it for spin pointers and assume unitary time evolution. In 
this case, the Hamiltonian of the coupling between the system and the n'th pointer 
measuring Bi is H n i = |Sj (gi 0^8 (t) (where e « 1) ; therefore the measurement 
propagator is 

U ni = e-^- B ^y ] = cos - i sin (|i? 4 ) ® ff W . ( 2 ) 

Assuming that both pointers initially point at \f x ) = ^= (|4-z) + |t*))j we denote the 
initial state as p in = J^Pa IV'a.sys, tx, tx> (V'a.sys, tx, tx|, where the label 'sys' stands for 

a 

'system'. The final state, p^, is given by 

pr = u (t) t/t (t) UyU(t) u liPin ul^ (t) U ] 2j U (t) C/t (t) 

= U (t) (cos {t Bj (t 2 )) - zsin {t Bi (t 2 )) g aM) 
x (cos (|s 4 (h)j -ism [^Bi ® crj, 11 ) p in 
x (cos (|Si (t x )) + zsin (|Si (ix)) ® <J l y ] ) 
x (cos (iBj (ta)) + *sin (1^ (t 2 )) ® a^) (t) , 

(3) 

where U (t) is the time propagator. The pointers' expectation value and variance are 
therefore 

E (crW ® 4 2 1) = ^Tr ({B z ( tl ) , ^ (f a )} p in ) + (e 4 ) , (4) 

^(4 1] ®4 21 ) =i + o( e 4 ) ■ (5) 

Similar calculations for continuous pointers are found at [BJ. 

It can be shown that following a single correlation measurement procedure (for 
short times), 

Pfi.sys = /Oin, S ys + C (e 3 ) , (6) 

i.e., the state is hardly influenced by the measurement, hence we refer to this 
measurement as "weak" [TJ [H] (Note that we do not use post selection). This is the 
key feature that enables this measurement method to work. As shown in equations 
([J) and ^ the "weakness" of the measurement comes with a price, every single weak 
measurement is highly inaccurate. This can be compensated with a large number of 
measurements. The precision of N such measurements of a single correlation scales 
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with N and e as Max | -3-^7^ , e 2 1 , hence the number of measurements required to 

reach a precision of S scales as N — O (S~ 4 ) . 

This procedure can be easily generalized to general (non-unitary) time evolutions 
and pointers of general dimension. 

3. Constructing the dynamics for Discrete Systems 

The Kraus representation of a superoperator is 

p(t) = J2KAt,to)p(to)Kl(t,t ) , (7) 

where the following constraint holds 

^Tt (Mo)jMMo) = 1 . (8) 

Eq. (|7|) determines the time evolution of density operators. Next we shall use 
this equation in order to derive the time evolution of the covariance matrix defined 
above. Then we show that by measuring covariance matrices one would estimate the 
dynamics. 

By transferring the time dependence from the Schrodinger picture to the 
Heisenberg picture, demanding that the expectation value of a general operator 
O remains unaffected by the picture, one can deduce the equivalent equation for 
operators: 

d(t)=J2 K U t ^)0(to)K fi (t,t ) . (9) 

Now, given a D level system, we wish to apply Eq. ([9]) to a complete basis of operators 
{B a } = {B ,Bi\ (i e {1..D 2 - l}), such that B = t D xD and {B. t } are chosen to be 
hermitian matrices. Since the set is complete the Kraus operators can be represented 
as 

K li (t,to)=u att {t,t )B a (10) 

(summation on double indices is assumed) where are complex vectors of coefficients 
and the following relations hold: 

[B a ,B b ] =if abc B c , (11) 
{B a ,B b } =g abc B c , (12) 

where f abc and g abc determine the structure of the basis. Substituting Eq. (fT0 |) -(|T2 | 
in Eq. ([9]) we obtain: 

B a (t) = ^U* bll U ca {if b ad + g b ad) {if dee + 9dce) B e (t ) , (13) 

where B a (t ) = B a . We set My(*,i ) = \ u tn u cii [ifud + 9ud) (ifdcj + 9dcj) and 
Xi{t,t ) = u% IJ u cll {ifm + 9bid)(ifdc0 + 9dco)- Substituting these definitions in Eq. 
(|13p we obtain 

B, (t) = My (t, t ) Bj {t ) + Xl (t, to) • (14) 
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Substituting Eq. (TT4)) in Eq. ([T]) we obtain 

crik (t, t ) = My (t, t ) a jk (t , to) ■ (15) 

This equation defines the time evolution of the covariance matrix. M could be 
retrieved from the covariance matrices simply by 

M(t,t ) = a{t,to)<J- 1 (t ,t ) . (16) 

In general a (to, to) might be singular. This happens if and only if the initial density 
matrix is singular (see appendix for a proof). In this case the density operator doesn't 
sample all the states in the Hilbert space, and so it is impossible to gain a complete 
knowledge of the evolution. However, it is worth mentioning that when picking a 
random initial state over a continuous uniform distribution, the chance of producing 
a singular density matrix is zero. 
From Eq. (fbT)) we get 

X (t, t ) = (B (t)) - a (t, t ) a' 1 (t , t ) (B (t )) . (17) 

Note that the estimation of \ (t, to) does not require measurements in addition to the 
ones needed in order to determine a (t, to). 

While M and x entirely encode the dynamics, it would be convenient to go back 
to the Kraus representation. In order to find a possible set of Kraus operators we use 
additional constraints on u aM 's that are derived by Substituting Eq. (fT0|) in Eq. ([8]) 

h u ln u bn(ifabc + 9abc)B c = l, (18) 

which means 

2 M aM U V (ifabk +9abk) = 0, (19) 

2 U *ay. U ^ (*/ab0 + 9abo) = 1 • (20) 

The symmetrical part of M, anti-symmetrical part of M and the trace of M, together 
with x an d the two constraints of Eq. (|T9j) and (j!Z0)) form 6 equations with 6 unknown 
objects which are u^u ^, Re (u^u 0fJ ), Im (u^uq^) , Re (u^Uj^) , Im (u*^u m ) and 
u *ifi u i^- The solution of this set of equations fully determines the matrix u M ® u'^ . 
Since this matrix is hermitian it can be diagonalized. Using the resulting sets 
of D orthonormal eigenvectors {v^} and corresponding eigenvalues {A M }, we get 
= y/X^v^ (no summation on /i), and thus a possible set of D 2 Kraus operators 
that govern the system's dynamics is given by 

= TVvSa • (21) 

Since this method relies on the measurement of o~ij (to, to) and (t, to) for every 
i,j G {l..D 2 — l} its complexity grows with D as C(_D 4 ). For n qudits (d level 
systems) D = d n and thus the complexity grows with the number of qudits as O (d 4 ™) , 
the same as in pQ. 

Notice that the described method works for general non-singular states, and in 
particular maximally mixed states and thermal states. An alternative procedure, 
based on covariance matrices defined at one time only, would not work with general 
states. 

Given the ability to measure M and x a t small time intervals it is possible to 
estimate the time derivative of the Kraus operators. For Markovian systems this 
allows one to estimate the Lindblad equation as described in [9]. 
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4. Constructing the dynamics for Gaussian channels 

In general, our proposed method requires the use of a complete set of operators. 
Nevertheless, in special cases it can be applied using a limited set of observables. One 
important example for this statement is the class of systems described by quadratic 
Hamiltonians H = Xijfjirjj + a^i (in this case 77 = (x,p) where x and p are arrays 
of coordinate and conjugate momenta respectively). The evolution of a subset of the 
system is then described as a Gaussian channel [lOiH"!"] that dictates for the conjugate 
coordinates an evolution of the type 

m (t) = My (t, t Q ) r)j (to) + X l (t, t Q ) . (22) 

This equation is identical to Eq. (|T4"]) . Therefore the rest of the derivation of the 
dynamics is the same as described in the previous section. (Note that in this case g a b c 
does not contribute to the calculations). While in the case of discrete systems the 
density matrix must be non singular in order for a (to, to) to be non singular, in the 
case of the Gaussian channels, the method would work for every state, even for states 
that do not sample the whole basis of states. The proof goes as follows: according to 
|12l 1 1 3| for every covariance matrix of canonical operators there exists a symplectic 
matrix S and a diagonal matrix W that satisfy: W = Sa (t,t) S T where Wa > \. 

Since S is symplectic it is invertible and so cr(t,t) — S~ l W (S T ) . In this form 

its clear that a (t^t)^ 1 — (S T ) 1 W~ 1 S always exists, in particular for t = to- The 
meaning of this feature is that in the case of Gaussian channels, our method is truly 
state independent. 

Since this method relies on the measurement of cry (to, to) and cry (t, to) for every 
i,j £ {1..2n} (where n is the number of particles) its complexity grows with the 
number of particles as O (n 2 ) . 

It is interesting to remark that in this Gaussian case, there is an analogy between 
the proposed method and Gaussian state tomography. In the "spatial problem" it 
is well known that the spatial correlations encoded in the covariance matrix fully 
determine a Gaussian state [13]. In our case we see that the Gaussian dynamics is 
encoded in the temporal correlations which form a temporal correlation matrix. 

5. Example: Constructing the Kraus operators for qubit dynamics 

For qubit dynamics we choose the natural operator basis {Bo = 1, -Bj = cr^}. The basis 
structure is therefore 

fijk = Ztijk , (23) 

9ab0 — 9a0b = 90ab = ^ab , (24) 

and the rest of the coefficients vanish. Using the explicit form of f a b c and g a b c we 
calculate 

My = (uq^Uo^ - U* k ^U kll ) Sij 

+ 21m (u* u k ) e ljk + 2Re (u*u ) , (25) 
Xi (t, to) = 2Re (u*oU t ) + Im (u*u k ) e ijk ■ (26) 
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Recalling that in this case D = 2 we obtain 

Tr(M) =2u* 0ll u 0ll , (27) 

Mij+Mji = 2(u* 0ll u 0t _ l -u* kll u kfl )S l:j +4:Re(u* ll u Jt _ l ) , (28) 

M l3 - Mji = 41m {u* 0p u kli ) e ijk . (29) 

From the constraints of Eq. (fTI))) . (f2T)|) we obtain two additional equations: 

2Re {u* Qli u kll ) - Im (u* lfl u ni ) e ijk = , (30) 

u 0n u 0n + u*^u ip , = 1 . (31) 

The solution of equations (f26)) - (|31[) is 

=^Tr(M), (32) 

Re (u^uon) = -Xi , (33) 

Im « M M o M ) = g e ijfc ( M fej - M jk ) , (34) 

Re = i(M ii +M ii + 2(1-Tr (M))5 y ) , (35) 

Im « M M J>) = l € ijkXk ■ (36) 



Equations (|3^|) - ([3^|) construct the matrix u M ® uj^ which is used to calculate a possible 
set of Kraus operators as explained above. This result can be easily generalized for n 
interacting qubits. 

Discussion 

In the present work we have proposed a method for implementing process tomography 
which is based on the information encoded in temporal correlations. We have shown 
that such correlations, embodied in the temporal correlations matrix, provide a general 
state-independent method to reconstruct the dynamical evolution in terms of Kraus 
operators. The complexity of our method grows exponentially with the number of 
particles in the discrete case (as in the usual approach), and quadratically with the 
number of particles in the Gaussian channel case. Since the proposed approach does 
not require the preparation and final state tomography, it is also applicable when the 
system is in a general (non singular) mixed state and in particular in a thermal state. 
In this respect, we believe that the present method is potentially useful for various 
experimental systems, such as for example in NEMS [TSl H2] and mesoscopic systems. 

We believe that in systems where the time scales are very short, and the coupling 
constant between the system and the measuring device is small it would be difficult 
to implement strong measurements. Since our method utilizes weak measurement, we 
believe our method would be suitable for such systems as well. Weak measurements 
have been recently realized in various systems [T71 [TB] . 

Finally, we note that it would be interesting to further investigate whether 
the present method can be combined with existing efficient methods for selective 
estimation of dynamics such as I3JS]. 
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Appendix: Proof that the covariance matrix is singular iff the state is 



Recall that the basis of operators is defined such that B — Idxd- First we prove 
that if the state p is singular then the covariance matrix a is singular: without loss of 
generality a singular state can be written in the form 



singular 




(37) 
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We choose {-Bj} such that 



Bi = 



Od-Ix-D-1 

1 



) 



(38) 



hence pB x = B 1 p = 0, and so 



o-u = Tr (p {B u Bi}) - 2Tr (pSO Tr (^BO 

= Tr (pBiBi) + Tr (SipBi) - 2Tr (pBi) Tr (pBj) 
= 0. 



(39) 



Since cij = an = for every i, a is singular. 

Now we prove that if a is singular then p is singular: if a is singular there must 
exist a certain Bi (which we refer to as B\) whose covariance with every member of 
{Bi} is zero, and in particular it's variance is zero. Therefore Bi| ran , ■, = c— number. 
Since B\ cannot be proportional to the unity range (p) < D, and thus p is singular. 



